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ABSTRACT: This paper proposes an uncertainty analysis framework based on the characterization of the un- 

certain parameter space. This characterization enables the identification of worst-case uncertainty combinations 
and the approximation of the failure and safe domains with a high level of accuracy. Because these approxi- 
mations are comprised of subsets of readily computable probability, they enable the calculation of arbitrarily 
tight upper and lower bounds to the failure probability. A Bernstein expansion approach is used to size hyper- 
rectangular subsets while a sum of squares programming approach is used to size quasi-ellipsoidal subsets. 
These methods are applicable to requirement functions whose functional dependency on the uncertainty is a 
known polynomial. Some of the most prominent features of the methodology are the substantial desensitization 
of the calculations from the uncertainty model assumed (i.e., the probability distribution describing the uncer- 
tainty) as well as the accommodation for changes in such a model with a practically insignificant amount of 
computational effort. 


1 INTRODUCTION 

This paper studies the reliability of a system for which 
a parametric mathematical model is available. The 
acceptability of the system depends upon its ability 
to satisfy several design requirements. These require- 
ments, which are represented by a set of inequality 
constraints on selected output metrics, depend on the 
uncertain parameter vector p. The system is deemed 
acceptable if all inequalities are satisfied. The require- 
ments/constraints partition the uncertain parameter 
space into two sets, the failure domain, where at least 
one of them is violated, and the safe domain, where all 
of them are satisfied. The reliability analysis of this 
system consists of assessing its ability to satisfy the 
requirements when the uncertain parameter p is free 
to take on any value from a prescribed set. The most 
common practice in reliability analysis is to assume a 
probabilistic Uncertainty Model of p (i.e., the random 
variable that models the uncertainty), and estimate the 
corresponding probability of failure. Calculating the 
failure probability is usually difficult since it requires 
evaluating a multi-dimensional integral over a com- 
plex integration domain. Sampling-based approaches 


(Niederreiter 1992, Kail and Wallace 1994) and meth- 
ods based on asymptotic approximations of the failure 
domain (Rackwitz 2001) are the engines of most (if 
not all) of the numerical tools used to estimate this 
probability. 

This paper proposes techniques that characterize 
the uncertain parameter space with a high level of fi- 
delity. A significant thrust of this research is the gen- 
eration of sequences of inner approximations to the 
safe and failure domains by subsets of readily com- 
putable probability. These sequences are chosen such 
that they almost surely fill up the region of interest. 
The strategies proposed, which are only applicable 
to requirement functions having an explicitly known 
polynomial dependency on the uncertainty, yield re- 
sults whose correctness is formally verifiable. The 
companion paper (Crespo et al. 2011) proposes strate- 
gies with the same goal but applicable to unrestricted 
requirement functions. Overall, the methodology en- 
ables the substantial desensitization of the calcula- 
tions from the assumed uncertainty model as well as 
the accommodation for changes in such a model with 
a practically insignificant amount of computational 
effort. 



This paper is organized as follows. Basic concepts 
are established in Section 2. Section 3 presents strate- 
gies for generating and refining approximations to 
the failure and safe domains. These approximations, 
along with the developments in (Crespo, Kenny, and 
Giesy 2011), enable the calculation of upper and 
lower bounds to the failure probability. Finally, a few 
concluding remarks close the paper. Proofs are omit- 
ted due space limitations. 

2 BASIC CONCEPTS AND NOTIONS 

Uncertainty models of p 6 I s can be probabilistic 
or non-probabilistic. A set whose members are all 
possible uncertain parameter realizations is a non- 
probabilistic model. This set, called the Support Set, 
will be denoted as A C M s . In a probabilistic uncer- 
tainty model, p is a random vector. This model is fully 
prescribed by the joint probability density function 
f p (p) : A — > M or the cumulative distribution func- 
tion F p (p) : A -»• [0,1]. 

Consider a system that depends on the uncertain 
parameter p. The design requirements imposed upon 
such a system are given by the vector 1 inequality 
g(p) < 0, where g : V — > !C, and ACPC ML The 
set D, where the constraint functions are defined, will 
be called the master domain. 

The failure domain, denoted as 'F C is com- 
prised of the parameter realizations that fail to satisfy 
at least one of the requirements. Specifically, the fail- 
ure domain is given by 

V 

^ = [j{P-9i(p) >°}- ( 1 ) 

2=1 

The safe domain, given by S — C(T), where G'( - ) de- 
notes the complement set operator given by C(X) = 
V\X, consists of the parameter realizations satisfy- 
ing all the design requirements. The failure probabil- 
ity associated with a probabilistic uncertainty model 
is given by 

P[?\ = Jf P (p)dp, (2) 

where P[-] is the probability operator. Techniques 
for approximating F and S will be presented below. 
The resulting approximations are comprised of hyper- 
rectangles or quasi-ellipsoids. 

The hyper-rectangle having m > 0 as the vector of 
half-lengths of the sides and p as its geometric center, 

'Throughout this paper, it is assumed that vector inequalities 
hold component-wise, super-indices denote a particular vector 
or set, and sub-indices refer to vector components; e.g., pj is the 
ith component of the vector p 1 . 


is given by 

1Z(P, m) = {p’.p — m < p < p + m} . (3) 

An alternative representation of this hyper-rectangle 
is given by 

H(p,m) = 5(p — m,p + m), (4) 

where 

5(x, y ) = [*!, yf\ x [x 2 , y 2 ] x • • • x [ x s , y s \, (5) 

is the Cartesian product of intervals. 

A subdivision is the process of dividing a set into 
subsets. Let p(-) be an operator whose input is any 
given set and its output are the subsets. A bisection- 
based subdivision in the ith direction is given by 

p(lZ) = {P(j p + w . m — w),1Z(p — w . rri — in)}, 

where w = [0, . . . , 0, m*/ 2, 0, . . . , 0]. Alternatively, 

p(ll) = {Afp 1 , v 1 + m), . . . , 5(v 2S ,v 2S + m)}, 

where v k is a vertex of 5(1,1 + m), leads to 2 s rect- 
angular subsets each of volume flL i rn i- 

The quasi-ellipsoid having m > 0 as the semi- 
principal axes vector and p as its geometric center, 
is given by 

£(p,m,n)= jp:(£(^) ) <lj <«> 

where n is an even natural number. Note that S is 
a closed set in R s having a polynomial boundary of 
degree n. Further notice that S ip. m, n) approaches 
7 Z(p, m) asymptotically from the inside as n oo. 

When applicable, p will be called the nominal pa- 
rameter point. The variable a will be used to deter- 
mine the containment of p in S. Specifically, a = 1 if 
g(p) < 0, otherwise a — — 1. 

The sections that follow provide the means to gen- 
erate and sequentially refine approximations of the 
failure and safe domains. 

3 REQUIREMENTS WITH KNOWN 
POLYNOMIAL DEPENDENCIES 

The key development in this section is the calcula- 
tion of inner and outer approximations to the fail- 
ure domain. These approximations are comprised of 
a collection of almost disjoint hyper-rectangles or 
quasi-ellipsoids. The almost disjoint condition is re- 
quired for estimating failure probability bounds (Cre- 
spo, Kenny, and Giesy 2011). 



Let T sub and S sub denote inner approximations 
(sub- sets) of the failure and safe domains. Thus, 

jrsup _ j s an ou ter approximation (super- 

set) of the failure domain. Note that 0 C T sub C T C 
jrsup q jy p ur th er notice that C(T sub U T sub ) con- 
tains the failure domain boundary dT . 

The strategies presented below can only be applied 
when the dependency of g on p assumes a known 
polynomial form. This dependency may occur natu- 
rally or may be the result of approximations. When 
the evaluation of g is computationally expensive, it is 
a common practice to use some realizations of the un- 
certain parameter and of the corresponding value of 
the constraint functions to build a surrogate model of 
g(p). If this model is chosen to be polynomial, the 
approaches introduced hereinafter can be deployed. 
Note however that the approximation error present in 
the surrogate model not only undermines the rigorous 
character of the methodology but may also yield re- 
sults that are invalidated by the actual g. 

Bernstein expansion-based techniques (Zettler and 
Garloff 1998) and interval analysis-based techniques 
(Jaulin et al. 2001) are well suited for hyper- 
rectangular geometries. Both of these approaches can 
be used to approximate sets with hyper-rectangles. 
The approximations resulting from the latter tech- 
nique however, are considerably more slack than 
those of the alternative technique for the same number 
of subsets. Sum of Squares (SOS) programing-based 
techniques (Packard et al. 2010) on the other hand, 
are better suited for sets with polynomial boundaries, 
e.g., ellipsoids. 

Formulations that enable evaluating set contain- 
ment, determining maximal deformations, and ap- 
proximating the failure and safe domains are con- 
sidered subsequently. By set containment we imply 
determining if all the members of a set are con- 
tained by another set. A maximal deformation (Cre- 
spo et al. 2008, Crespo et al. 2009) is a process 
that enables quantifying the separation between a 
point and dT as well as the identification of worst- 
case uncertainty combinations. By failure or safe do- 
main approximations we imply the generation of se- 
quences of inner and outer approximations to such 
sets. The sequence jF| ub , . . .} contains inner 

approximations of the failure domain. The sequence 
{T{ up ,T S ^ V , ...} contains outer approximations of 
the failure domain. Similarly, . . .} and 

{C(Tl ub ), C(T £ ub ), . . .} are sequences of inner and 
outer approximations to the safe domain. All these 
sequences approach the domain being approximated 
when their number of terms increases. 


4 BERNSTEIN EXPANSION 

The Bernstein expansion approach (Zettler and 
Garloff 1998) enables determining if a polynomial 
inequality holds or not over a hyper-rectangular do- 
main. This approach requires mapping this domain to 
the unit hyper cube. Denote by u = U (p) an affine 
transformation that maps the hyper-rectangle V onto 
the unit cube U = 7 £(|, |). Then h(u) = g(l 7 _1 (w)) 
is a polynomial on U. 

For simplicity in the presentation we first consider 
the case where there is a single constraint function 
assuming a univariate polynomial form. If the poly- 
nomial h is given by 

71 

h(u) = (7) 

i = 0 

its Bernstein expansion is 

n 

h(u) ='£b i (V,g)B?(u), (8) 

i = 0 

where 

B,”M = <9) 

is the zth Bernstein polynomial of degree n (i.e., an 

element of the basis) and 

i 

f>.(©.s) = £4A, do) 

3=0 G' 

is the ith Bernstein coefficient. Some fundamental 
properties of this basis are B™(u) = 1 (normal- 
ization), 0 < B™(u) < 1 for 0 < u < 1 (bounded- 
ness), and B)'(u) = ( 1 — u) > 0 (symmetry). 

Simple manipulations lead to the free function eval- 
uation property: 

h(0) = b 0 (V,g), (11) 

h(l) = b n (V,g). (12) 

The range enclosing property follows directly from 
the normalization and boundedness properties and is 
given by 

min {bi{T>,g)} < g(p) < ma x.{bi{V,g)} , (13) 

i<n i<n 

for all p G D. Note that the Bernstein expansion en- 
ables bounding the range of the polynomial by mere 
algebraic manipulations. In contrast to the nonlinear 
optimization approach (Crespo et al. 2011) and SOS 



approaches, there is no need to solve an optimization 
problem or even evaluate the polynomial. 

The single constraint, multivariate polynomial case 
is considered next. Define the multi-index i as a vec- 
tor of non-negative integers i — [ii, . . . ,i s ]. Note that 
monomials can be represented as u 1 = u] 1 u l f ■ ■ ■ u l s s . 
An .s-variatc polynomial can be represented as 

h(u) = y^djU 1 , (14) 

i<n 

where u eW. Recall that vector inequalities should 
be interpreted component-wise. The Bernstein expan- 
sion of (14) is given by 

h{u) = Y J bi{V 1 g)B?{u) } (15) 

i<n 

where 

B?(u) = BZ( Ul )---B?;(u s ). (16) 

is the ith Bernstein polynomial of degree n and 

n F. k \ 

b i (V,g) = ^2H-^-a j , (17) 

j<ik = 1 (/) 

«/ k 

for i < n, is the 7th Bernstein coefficient. The normal- 
ization, boundedness and symmetry properties extend 
to the multi-variate case. The free function evaluation 
property becomes 

h ([ii/rii, . . . ,i s /n s ]) = bi(V,g), (18) 

where i is an element of (0, r^} x • • • x (0, n s }. The 
range enclosing property becomes 

mm{bi(V,g)} < g{p) < max{ bi{V,g)}, (19) 

i<n i<n 

for all peD. Tighter bounds on the range of g(p) 
are obtained by subdividing V into subsets, calculat- 
ing the Bernstein expansion for each of them, and se- 
lecting the largest and the smallest of all Bernstein 
coefficients. In particular, if p(D) = {TZi , . . . ,lZ t } is 
a subdivision of the master domain, then for all p EV, 

min < min {bi(R, r g)} > < g < max < max {bi('R r g)} 

7 ^ t I z R I 7 ^ t I 1^71 

( 20 ) 

These bounds converge to the global minimum and 
global maximum of g(p) when the volume of the sub- 
sets in p approaches zero. The reason the bounds in 
(20) can be calculated efficiently is that the Bernstein 
coefficients on the elements of p(T>) resulting from 
applying the bisection-based subdivision of Section 2 
can be calculated directly from those on D (Zettler 
and Garloff 1998). The sequential application of the 
bisection-based subdivision scheme leads to bounds 
that converge to the global minimum and maximum. 


4. 1 Set Containment of Hyper-Rectangles 

The following theorem uses the developments of the 
previous section to determine whether a set TZ(p, m) 
is fully contained in the safe domain S or failure do- 
main F . 

Theorem 1. Given a subdivision p( V) = 
(77i, . . . the set containment condition 

lZ(p, m) C S holds if 

max < max < max {bfIZj , g fe )} > > < 0. (21) 

k<v I j<t I i<n I J 


Furthermore, 'Rip, m) % S if there exists k < v, j < 
t, and a multi-index i E (0, rii} x • • • x (0, n s } such 
that 

bifRj,g k )> o. (22) 

While Formula (21) results from choosing the 
largest upper bound in (20) over all the constraint 
functions, Formula (22) results from applying the free 
function evaluation property. 

Theorem 2. Given a subdivision p(P) = 

{TZi, . . . the set containment condition 

TZ(p, m) C F holds if 

max < min < min {bi(TZj,g k )} > i > 0. (23) 

k<v I j<t I i<n \ \ 


Furthermore, Rip, mi) f F if there exists j < t, and 
a multi-index i { 0 . rr i } x • • • x { 0 . n s } such that for 
all k < v, 

bFR, r g k ) < ()• (24) 

While Formula (23) results from choosing the 
largest lower bound in (20) over all the constraint 
functions. Formula (24) results from applying the free 
function evaluation property to any of them. 

The asymmetry between Formulas (21) and (23) re- 
) suits from the definitions of the safe and failure do- 
r mains. In order to apply these set containment condi- 
' tions one is required to implement an algorithm for 
sequentially subdividing the master domain. When 
checking whether 7 Z(p,m) C S holds, the subdivi- 
sion algorithm should be stopped when either (21) or 
(22) is satisfied. When checking whether 7 Z(p, m) C 
F holds, the subdivision algorithm should be stopped 
when (23) or (24) is satisfied. The strategies in (Smith 
2009) used to compute Bernstein coefficients were 
adopted. The computational complexity of the re- 
sulting subdivision algorithm is nearly linear with 
the number of monomials. This makes the Bernstein 
polynomial approach very efficient. 



4.2 Maximal Deformation of Hyper- Rectangles 

In this section we use the Bernstein expansion ap- 
proach to evaluate the set containment conditions re- 
quired to perform homothetic deformations. In partic- 
ular, one can solve for 

a = sup {a : 7 Z(p,am) C 5}, (25) 


4. Let P*^f = C(S-ff). Evaluate P[Pff\] and 
P[P-^\ (Crespo et al. 2011). 

5. If P[Pff b ] > 1 — P m ax declare the system ac- 
ceptable and stop. If P[P^\\ < P m .ax declare 
the system unacceptable and stop. Otherwise, in- 
crease i by one, and go to Step (1). 


by evaluating the set containment condition with For- 
mula (21) after replacing 7 Z(p,m) with 7Z(p, am). 
An identical process can be applied to compute 

a = sup {a : 7 Z(p, am) C JF}, (26) 

using Formula (23). 

Alternatively, one may search for the maximal de- 
formation using nonlinear optimization (Crespo et al. 
2011) and then use the developments in Section 4. 1 to 
verify convergence to the global optimum, e.g., show 
that 7 Z(p, am) C S. 

4.3 Failure Domain Approximations 

The algorithm below iteratively generates indexed 
sets A j, S- ub , and P? ub of hyper-rectangles where Sf uh 
is an inner approximation to the safe domain, P? ub is 
an inner approximation to the failure domain, and A* 
is a region whose containment in P or S is to be de- 
termined. The application of the set containment con- 
ditions (21) and (23) enable the expansion of the inner 
approximations. The algorithm proceeds by selecting 
a hyper-rectangle 1Z e A, : . If 1Z C S, the inner approx- 
imation to the safe domain is expanded with this ele- 
ment. If 7Z C P , the inner approximation to the failure 
domain is expanded with this element. If none of these 
conditions are satisfied, the elements of a subdivision 
of this rectangle are appended to A*. The algorithm 
terminates when the bounds to the failure probability 
exceeds a prescribed limit. The algorithmic represen- 
tation of this procedure is as follows. 

Fet g(p) < 0 denote the set of system requirements 
and f p (p), p G A C V be the uncertainty model. De- 
note by Pmax a parameter prescribed by the user de- 
noting the largest admissible failure probability asso- 
ciated with the system. Set i — 1, Ai = {72}, P( ub = 0, 
and S{ ub — 0. 

1. Fet 7Z be a largest element of A*. 

2. Perform a Bernstein expansion of all constraint 
functions in 7 Z. 

3. If Equation (21) holds let A i+ i = A* \ 7 Z, \ = 

S[ ub U 7 Z, and P t f b = Pf ub . If Equation (23) 
holds let Ai+i = A ?; \ 77, Pff\ = P t sub U 77, 

\ = S[" h . If neither equation holds, set A ?;+1 = 
(A * \ FZ) U p(1Z), Sff b = S° ub , and Pff b = P s t ub ■ 


As the number of iterations increases, Sf ub and 
Pf" h approach the safe and failure domain. P\P:f ,:b ] 
and P[S? ub ] are monotonically increasing functions 
of i. On the other hand, P[Aj] and P\P slip ] are 
monotonically decreasing functions of i. Note that 
the elements in A, : are an approximation of OP. The 
larger the value of i the better the approximation. 


Example 1: Consider the constraint functions 


9 1 


2 4, 4 2 o 2 2 

P 1 P 2 + P 1 P 2 ~ 3PiP 2 ~ P 1 P 2 + 


p\ + pI 
200 


7 

loo’ 

(27) 


92 = -p\p\ - p\p\ + 3Pi/?2 + - °- 9 > (28) 

These requirement functions have been chosen so the 
safe domain is multiply connected. Figure 1 shows 
Pf ub and S- ub for a fixed value of i. In this case we 
have used the bisection-based subdivision where the 
variable being subdivided alternates. The set of boxes 
for which set containment cannot be established are 
colored in white. Note that this region is a tight ap- 
proximation of OP. Tight bounds to the failure prob- 
ability can be readily calculated from these approxi- 
mations. 


5 SOS PROGRAMMING 

Even though the theory supporting the developments 
that follow applies to all sets having polynomial 
boundaries, we focus on quasi-ellipsoids since they 
enable the analytical calculation of their probability 
(Crespo et al. 2011). 

5 . 1 SOS Determination 

For simplicity in the presentation we first consider the 
single constraint case. If g(p) is a polynomial of de- 
gree less than or equal to 2d in the variable p e M ,s , its 
vectorial representation is 

g(p) = c T x(p), (29) 

where x : — > M dim ( a; ) is a vector of monomials in 

p of degree less than or equal to 2 d, and c is a vector 




hi,...,h n such that g = V'' , If. The polynomial 
g(p) is a SOS if and only if there exists a posi- 
tive semi-definite matrix Q, to be denoted as 0 A 
0, that satisfies (30) (Powers and Wormann 1998, 
Parrillo 2000). The functions hi,..., h n , that consti- 
tute the SOS representation of g, result from mak- 
ing Choleski or Schur decompositions of Q. Conse- 
quently, g is a SOS if an only if there exist a A for 
which Q o + Y^=\ A 0. This is a Linear Matrix 
Inequality (LMI) feasibility problem. 

Numerical techniques for solving semi-definite 
programs can be used to find a solution to this LMI 
feasibility problem. Publicly available software, such 
as SOSTOOLS, YALMIP and SeDuMi, automate the 
process of posing and solving this convex optimiza- 
tion problem. Unfortunately, its computational re- 
quirements (e.g., the number of monomials that re- 
quire representation and the dimension of the null 
space) grow rapidly with the degree of the polyno- 
mials and the dimension of the parameter space s. 


Figure 1: Tf ub (red), Sf ub (green), dT (thick line) and individ- 
ual constraint boundaries (thin lines). 


of coefficients. On the other hand, the Gram matrix 
representation of this polynomial is 

g(p) = z T (p)Qz(p), (30) 

where z : —> M dim U) is a vector of monomials 

in p of degree less than or equal to d and 0 e 
^dim(z)xdim(z) j s a S y mme tric matrix. The Gram rep- 
resentation of g is not unique. A procedure to param- 
eterize all possible representations is as follows. De- 
fine the linear operator L that maps each symmetric 
matrix Q to the polynomial coefficients c, i.e., 

L(Q) = c. (31) 

A matrix representation of L can be computed since 
both its domain and its range are finite dimensional. 
This transformation enables us to parameterize the 
family of symmetric matrices yielding Gram repre- 
sentations via 

p 

Q = Qo + J2 x * Ni ’ (32) 

i= 1 

where Qo is the symmetric matrix corresponding to a 
particular Gram representation (i.e., L(Q 0 ) = c), the 
set { Ni , . . . N p } is a basis of the null space of L (i.e., 
L(Ni) — 0 for i — 1, . . . ,p) and A G W is a vector of 
multipliers. Note that for any value of A in (32), (30) 
is a valid Gram representation of g. 

The Gram representation of a polynomial enables 
us to determine if a polynomial is a SOS. The poly- 
nomial g(p) is a SOS if there exist polynomials 


5 .2 Set Containment of Quasi-Ellipsoids 

We want to determine if the set £{p,m,n ) is fully 
contained in the safe or failure domains by using the 
developments of the previous section. These set con- 
tainment conditions, which are equivalent to 

£(p,m,n)QS, (33) 

£(p, m,n) C T, (34) 

can be evaluated using the following Theorem. 

Theorem 3. Let £ = {p : e(p) < 0}, where e(p) 
is polynomial, be an alternate representation of the 
quasi-ellipsoid in (6), I = (1, . . . , v} and X = {i £ I : 
g, (p) >0}. If there exist a function q(p) : M ,s — > K ,; 
such that q(p) > 0 and q i (p)e(p) — gfp ) > 0 for 
all i e I and all p G A then £(p, m , n) C S. On the 
other hand, if there exists a function q(p) : M s — > M 
such that q{jp) > 0 and q{p)e(p) + gfp) > 0 for some 
i G X and cdl p € A then £ (p, m , n) C T. 

Note that there are 2 positive semi-definite (PSD) 
constraints per constraint function. PSD constraints, 
even polynomial ones, are numerically difficult to 
solve. However, by restricting the constraint func- 
tions g to take on a multivariable polynomial form; 
and relaxing the constraints to be SOS polynomials, 
in which case Theorem 3 still holds, we obtain a 
problem that can be solved efficiently. While all SOS 
polynomials are PSD, not all the PSD polynomials 
are SOS. Therefore, the relaxation to SOS polynomi- 
als yields a sufficient condition for the global non- 
negativity constraints in Theorem 3. 

In practice, the SOS polynomial multipliers q(p) 
and q(p) must be restricted to be in a fixed finite 




dimensional subspace of polynomials (e.g., quartic 
polynomials). This is achieved by prescribing a poly- 
nomial basis for them, and using semi-definite pro- 
graming to search for the coefficients of the cor- 
responding linear combination. Note however that 
the prescription of an overly restrictive basis may 
wrongly invalidate the set containment condition 
(e.g., a quadratic basis for q is prescribed but a quar- 
tic q is needed to demonstrate containment). Fur- 
ther notice that the SOS approach cannot be used to 
demonstrate that there is no set containment. Numer- 
ical experiments indicate that problems on the order 
of s + n + deg(< 7 j) < 12 are solvable. 

5.3 Maximal Deformation of Quasi-Ellipsoids 

In this section we use the SOS programing approach 
to evaluate the set containment conditions required by 
the homothetic deformations. In particular, we want to 
solve for 


a = sup{« : £(p,am,n ) 

C5}, 

(35) 

a = sup{a : £(p,am,n) 

cx}. 

(36) 


Theorem 3 enables the reformulation of (35) as 
a — ( max 8 | subject to (37) 

VT <J(p) J 

q(p) G 9[p\, (38) 

~9i(p) + Qiip) (e(p) -P + 1) e. 9[p\, Vi G /, (39) 

a — p < x/ \J~Pm < b — p, (40) 

where A = S(a,b), and 9\p] denotes the set of SOS 
polynomials in the variable p. Constraint (40) ensures 
that the maximal set is fully contained in the support 
set so its probability can be calculated analytically. 
The application of Theorem 3 to (36) yields 

p = max p t , (41) 

iex 

where 

Pi — [ max P ) subject to (42) 

\0, g{p) J 

q(p) G 9\p\, (43) 

9 tip) + q(p) Hp) - P + !] e 9[p\, ( 44 ) 

and constraint (40). In general p < a. This is so be- 
cause £ (p, pm, n) will be limited by one of the man- 
ifolds gfp) = 0 with i G X, but such a manifold may 
not be part of dX. Therefore, while the formulation 
for containment by the safe domain always converges 


to the critical similitude ratio a, the one for contain- 
ment by the failure domain may only underestimate 
it. 

The above problems, which have an objective func- 
tion that is linear in the decision variables (i.e., 8) 
as well as multiple SOS constraints, are called SOS 
programs. The terms q,P in (39) and qP in (44) 
make both SOS programs bi-linear in the decision 
variables. The constraints containing these terms are 
quasi-convex; e.g., for a fixed value of P all con- 
straints are convex in q. Hence these problems can 
be solved sequentially by searching for max(P) via 
bisection and by searching for q via linear SOS pro- 
graming techniques. 

With the critical similitude ratio (or its estimate) in 
hand, the maximal set (or its estimate) can be read- 
ily calculated. The corresponding critical parameter 
value p (or its estimate) however, cannot be calculated 
directly from the solution to the above problems. 

5.4 Failure Domain Approximations 

The algorithm proposed here is similar to the one in 
Section 4.3 since it iteratively generates the indexed 
sets Sf ub (inner approximation to the safe domain), 
Tf ub (inner approximation to the failure domain), and 
Aj (region whose containment in T or S is to be de- 
termined). At any given iteration we first chose a rect- 
angle from A j. By the means presented in Section 5.2 
we determine if the quasi-ellipsoid inscribed in this 
hyper-rectangle is contained in the safe or failure do- 
mains. If the ellipsoid is contained in the safe domain, 
the inner approximation to the safe domain is ex- 
panded with this element. If the ellipsoid is contained 
in the failure domain, the inner approximation to the 
failure domain is expanded with this element. Other- 
wise, the rectangle is subdivided into smaller rectan- 
gles (see section 2 for two subdivision logics), and 
these subsets are annexed to A*. The iteration of this 
procedure, which is stopped when the bounds to the 
failure probability exceeds a prescribed limit, leads to 
the desired sequence of approximations. The algorith- 
mic representation of this procedure is as follows. 

Start from the very same initial iteration values of 
the Algorithm in Section 4.3. 

1 . Let 7 Z(p, m) be a largest element of A*. Set H = 

£(p,m,n). 

2. Solve Equation (37) if a = 1. Otherwise solve 
Equation (42). Let r denote either a or p. 

3. If r < 1, set Aj +1 = (A* \ K) U p(K), = 

S- ub , and JEf+i = lFf" b ■ If r > 1 and a = 1 let 
Aj +1 = A* \ 7^, Sff\ = Sf ub U n and = 
Tf lh . If r > 1 and a — —l let A i+ i = Aj \ 77., 

JFf+i = Tt ub U Q and S™ b = S, sub . 
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Figure 2: F- ub (red), Sf ub (green), dT (thick line) and individ- 
ual constraint boundaries (thin lines). 


4. Let Evaluate P[T^\\ and 

P\^i± i\ or their lower bounds and 

V;(^D depending upon the applicable Theo- 
rem (Crespo et al. 201 1). 

5. If P\P™\] > 1 — Pmax declare the system ac- 
ceptable and stop. If P[Ft+\] C Pmax declare 
the system unacceptable and stop. Otherwise, in- 
crease i by one, and go to Step (1). 


The closing remarks of Section 4.3 apply here as well. 

Example 2: The very same problem setup of 
Example 1 is considered here. Figure 2 shows the 
subsets comprising T sub and S sub for a fixed value 
of i. Note that the approximations, which result 
from uniting quartic polynomials (n = 4), have 
voids among neighboring subsets. Further notice 
that subsets of J r £"\ such as those in the vicinity of 
p = [-1,-1. 25] , could have been replaced by a 
single subset that contains them. Even though this 
larger subset is contained in T, the containment 
condition (43-44) is not satisfied. This is the mani- 
festation of ft < a. Bounds to the failure probability 
can be readily calculated from the approximations. 
Due to the empty space among the neighboring sub- 
sets comprising the approximations, these bounds, 
whose slackness decreases as n increases, are not 
as tight as those for approximations comprised of 
hyper-rectangles. 


This paper presents an uncertainty analysis frame- 
work applicable to systems subject to polynomial re- 
quirements. Approaches based on Bernstein expan- 
sions and SOS programing are proposed. These and 
all other methods requiring the exploration of the 
uncertain parameter space suffer from the curse of 
dimensionality, hence their computational demands 
grow exponentially with the number of uncertain pa- 
rameters. Unfortunately only this space provides the 
sense of causality required to understand and prevent 
failure. The high dimensionality of this space along 
with the inability to guarantee that optimization prob- 
lems posed there will converge to the global optimum 
are the main liability of the engineering decisions sup- 
ported by the outcomes of these methods. This paper 
proposes techniques that eliminate this liability when 
the dimension of the uncertain space is moderate. The 
algorithms proposed allow for data parallelism (i.e., 
perform computations simultaneously on elements of 
a subdivision of the master domain). This will help 
mitigate the formidable challenges of having a large 
number of uncertain parameters. 
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